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ABSTRACT 



We reconsider the collapse of cosmic structures in an Einstein-de Sitter Universe, using 
the self similar initial conditions of Fillmore & Goldreich (1984). We first derive a new 
CO , approximation to describe the dark matter dynamics in spherical geometry, that we refer to the 

"fluid approach" . This method enables us to recover the self-similarity solutions of Fillmore 
^ ■ & Goldreich for dark matter. We derive also new self-similarity solutions for the gas. We 

, thus compare directly gas and dark matter dynamics, focusing on the differences due to their 

I different dimensionalities in velocity space. This work may have interesting consequences for gas 

■ and dark matter distributions in large galaxy clusters, allowing to explain why the total mass 

, profile is always steeper than the X-ray gas profile. We discuss also the shape of the dark matter 

density profile found in N-body simulations in terms of a change of dimensionality in the dark 
matter velocity space. The stable clustering hypothesis has been finally considered in the light 
of this analytical approach. 



Subject headings: Cosmology: theory - dark matter - hydrodynamics 



INTRODUCTION 



X 

' The study of cosmic structures formation within the gravitational instability picture has proven to be 

a usefuU tool over the last decade to connect the initial power spectrum with the density profiles in dark 
matter halos. Gunn & Gott (1972) derived analytical density profiles obtained during the collapse of a 
self-similar dark matter density perturbation. This so-called secondary infall theory was further extended 
and refined by Fillmore & Goldreich (1984) and Bertschinger (1985). The latter also included a detailed 
treatment of gas collapse, and showed that gas and dark matter distributions are similar in the special case 
that he studied. Hoffman & Shaham (1985) generalized these works to Gaussian random fields, showing 
that for scale-free power spectra {P{k) oc fc"), dark matter halo are well approximated by a singular 
isothermal sphere (p oc r^^) for —3 < n < — 1, and that the density profile steepens for larger value of n. 
N-body experiments confirmed roughly these conclusions (Frenk et al. 1985; Quinn, Salmon & Zurek 1986), 
by taking precisely into account the complex dynamical behavior of 3D hierarchical clustering. 

Observational estimates of gas and dark matter mass profiles in X-ray clusters have improved over 
the past years by the emergence of new technics like strong and/or weak lensing analysis and temperature 
maps obtained by the X-ray emission. Two challenging properties of the observed density profiles in X-ray 
clusters are first that there is a significant segregation between gas and dark matter (David, Jones & 
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Forman 1995), and second that mass estimates in the inner part of clusters are not fitted by a singular 
isothermal sphere, but rather by a much shallower profile, like a r^^ law (Wu & Hammer 1993) or a very 
small core radius (Soucail & MeUier 1994). Very recently, Navarro, Frenk & White (1996) showed that it 
was indeed the case in high resolution N-body simulations (see also Pfitzner 1996). They obtained a very 
general density profile well fitted by the formula 



'^'^ " r(i 

with a single shape parameter r^. 

For the gaseous component, the situation is rather unclear. Numerical calculations performed by 
Pearce, Thomas & Couchman (1994) revealed that a small segregation between gas and dark matter occurs 
within the hierarchical clustering scenario. This would be due to a systematic energy transfer between gas 
and dark matter during mergers of the small sub-clumps that lead to the final cluster. On the other hand, 
Anninos & Norman (1996), by using a totally different numerical method, have studied gas and dark matter 
density profiles of a Coma-like cluster in a CDM cosmogony. They found no significant segregation between 
gas and dark matter. Moreover, they found that the dark matter density profile is not fitted by a Navarro 
et al. (1996) model, but rather by a unique power law p oc r~^/^. 

In this paper, we address these questions using an analytical approach. This has been possible only 
within a spherical geometry approximation, together with self-similar initial conditions. However, we are 
aware that spherical collapse is a strong approximation, because it suffers for example from the so-called 
radial orbit instability (Hcnon 1973; see also Binney & Tremaine 1987). Our work mainly extends the 
work of Fillmore & Goldreich (1984) and Bertschinger (1985) to the case of gas collapse, and discuss 
the consequences to more realistic initial conditions. We briefly recall the initial conditions and the 
basic notations used in Fillmore & Goldreich (1984) to derive their self-similar solutions. As we restrict 
ourselves to the spherical case, we simplify the original notations of Fillmore & Goldreich, by taking their 
dimensionality n = 3. 

Each gas or dark matter shell is labeled by the mass Mj initially enclosed within its initial radius at 
the initial time t,. The initial velocity field is supposed to follow the Hubble law 



(1) 



The initial perturbation in mass is assumed to be a power law of the initial mass 



Mi [MoJ 

Each spherical shell is expanding until the turn-around epoch, and the trajectory is described by 
Kepler's law until it crosses another trajectory. The turn-around radius and the turn-around epoch of a 
given shell Mj are given by 



r* 



e-l-l/3 



(2) 
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Another useful! quantity is the current turn-around mass which labels the current turn-around shell 

Fillmore & Goldreich (1984) derived self-similar solutions for the gravitational collapse of dark matter, 
using these rather peculiar initial conditions. Bertschinger (1985) computed for the particular case e = 1 
the solution for both gas and dark matter. In this paper, we derive both gas and dark matter self-similar 
solutions for < e < 1. In order to compare theoretically the dynamics of each component, we use a fluid 
approximation to describe the dark matter collapse. We show in the first section that this approximation 
is valid, because, as stated by Fillmore & Goldreich (1984), the mass M{r,t) is an adiabatic invariant. 
Using this fluid approach, we derive in section 2 the well-known Fillmore & Goldreich (1984) self-similar 
solutions in a totally different way than the original derivation. This gives us new insights in the nature 
of dark matter in spherical geometry, treated here as a one-dimensional fluid. We derive also in section 2 
self-similar solutions for the gas, and we find that they differ from the dark matter solutions. In section 
3, we present numerical simulations to validate our analytical treatment of gas dynamics, and finally, we 
discuss possible applications of our work to other observational or numerical studies. 



2. A FLUID APPROACH FOR DARK MATTER 

In the original paper of Fillmore & Goldreich, the relaxed region of the forming object, called the halo, 
has a mass profile parameterized by 



M{r,t) ^ nity (4) 
where the function of time k was parameterized by 



K{t) oc (5) 

This is justified by the scale free nature of the problem, and this is likely to be valid only in the most 
inner part of the halo, where the mass profile has reached asymptotically a power law both in radius and in 
time. Dark matter particles in the halo are then following the equation of motion 

^ = -— = -G-ity' (6) 

Fillmore & Goldreich (1984) main hypothesis was to assume that for particles deep inside the halo, 
the function K{t) was a slowly varying function of time compared to one orbit period. This means that the 
right-hand side of equation (^ can be taken as a function of radius only over one orbit period, and thus we 
can integrate equation (|[) to obtain an energy integral along each trajectory 
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d-t) ) 

where ra{t) is the apapsis radius of the current orbit. As we show it now, such an energy integral aUows to 
close the Vlasov hierarchy which describes the collision-less dynamics of dark matter. 

We define here the distribution function / in phase space of the dark matter particles, considering only 
radial orbits. Therefore, the distribution function / is defined as the mass of dark matter particles per unit 
radius r and per unit ID velocity u, and satisfies the Vlasov equation 



df df GMdf „ 
ot or ov 

The fluid limit of such a kinetic description is obtained by taking the n*'' moments of equation (|^) 
in velocity space. We then obtain for the first three orders the mass conservation equation, or continuity 
equation 



t + |:(,») = o (<,) 



the momentum conservation equation, or Euler equation 



dt ' ' dr ^ ' ' r2 
and the energy conservation equation 



ifiu) + 7^ (n + fiu^) = -M— (10) 



The mass per unit length is noted /i, the mean velocity of dark matter particles in a given fluid element 
< u > is noted u, and 11 is related to the velocity dispersion of dark matter particles in this fluid element 
by n = fi{< > —V?). Therefore, 11 is the analogous of a thermal pressure in the dark matter fluid. This 
hierarchy is a priori not closed, and one has to know the exact form of the third moment 6*3 to solve the 
problem. 

Using the energy integral (eq. 0), it appears that the fluid element at position r is crossed by particles 
with opposite velocities which follow trajectories labeled by their apapsis. The distribution function is 
therefore even in velocity space, and moments of odd order vanish. This implies that S'3 and also u ~ 0. 

The fluid approximation (5*3 = 0) allows us to close the set of equations (p[]ll|) which becomes the 
fluid equations for dark matter. Since these equations are hyperbolic, discontinuities in the flow appear 
naturally. They are treated using the Rankine-Hugoniot discontinuity relations. Physically, bulk kinetic 
energy l/2^v? is thus dissipated into internal kinetic energy 1/2/i < >. Relaxation in the dark matter 
fluid is equivalent to shell crossings and the discontinuity is analogous to the first caustic in the dark 
matter density profile. In the exact treatment of Fillmore & Goldreich (1984), relaxation occurs in a thick 
layer, between the first caustic (where S'3 differs significantly from zero) and the radius where the fluid 
approximation is roughly recovered (S'3 ~ 0). 



3. GAS AND DARK MATTER SIMILARITY SOLUTIONS 



Our fluid approach gives new insight in dark matter dynamics. It also allows us to compare directly 
gas and dark matter dynamics. Because in our fluid approximation we don't track particles trajectories, it 
is possible to derive the self-similar solutions in a rather simple analytical way. In fact, the mass enclosed 
by a given Lagrangian fluid element remains constant in time. Each fluid element is initially labeled by the 
initial enclosed mass Mi. Because this mass remains constant, we drop here the subscript i and we label 
each fluid element by its enclosed mass M . 

We first consider the case of the pure dark matter collapse, and recover Fillmore & Goldreich results, 
and finally we investigate the pure gas collapse, showing that the halo profile in both cases can differ 
significantly. 

3.1. Dark Matter Halo Profiles in the Fluid Approximation 

Each fiuid element first expands, and reaches turn-around at t = After turn-around, it falls towards 
the center and reaches the relaxation front. Self-similarity implies that this front is located at a constant 
fraction of the turn-around radius 

Tg — As 7*^ 

and that the fluid element get shocked at a constant fraction of the turn-around epoch 

Behind the shock, the fluid element gradually reaches the asymptotic regime. During this self-similar 
regime, we parameterize its evolution by 

r{M,t) (X MPt-1 (12) 

where q > 0. To ensure self-consistency, we extrapolate this regime back to t = tg. We then obtain using 
equations ^ and (||) 

13 

P=3+e+29e (13) 

Mass conservation implies that 7 = 1/p (eq. [Q) and s = —q/p (eq. |^). We then assume that each 
fluid element is in hydrostatic equilibrium. In the previous section, we noticed that this assumption was 
valid during the asymptotic halo regime. Note however that it doesn't imply necessarily that q — 0. We 
will see in the foUowings, that it is in fact possible to obtain a quasi-hydrostatic flow, which is not strictly 
stationary. 

Hydrostatic equilibrium writes then for a fluid element M 
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Integrating this equation between M and the currently shocked fluid element Ms(t), and using the scaled 
variable x — M/Ms{t), leads after a few algebra to 



n(M,t) = ns{t) + /' x'-'^dx (15) 



Rsit) 



M/M, 



which is consistent with the Virial theorem. The integral in equation ( [T5| ) converges or diverges whether 
p < 1 or p > 1. 

Let us assume first that p < 1. In the asymptotic regime, one has M/Ms{t) <C 1. The integral in 
equation (|l^) is thus a constant to leading order. Consequently, in the halo 11 is uniform, and is determined 
entirely by the post-shock conditions. The hydrostatic equilibrium equation cannot be satisfied in that case. 
The validity range for p is therefore p > 1. We deduce then the validity range for q from equation (p^. 

e<? q>lC^T;-e) (16) 



3 " - 3e V3 

As p > 1, the integral diverges, and to leading order, one gets 



RS? \Ms{t) 

It is straightforward to check that the hydrostatic condition (eq. |l4j) is fulfilled in that case. We need 
another equation to find the final solution of the problem. The last equation on 11 states that twice the 
internal energy of a given shell is equal to its gravitational energy. The action is thus given by 



Integrating by part and using again the hydrostatic equilibrium assumption leads to 



t ^ 



Since > 0, the action is always a constant in the asymptotic regime where t S> i*, and takes the value 

S{M) = q'l (17) 

Applying the least action principle means here that the solution is obtained for the minimum value of q. 
From equation (p^, we find 

2 
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In this way, we finally recover the solution of Fillmore & Goldreich within our fluid approximation. The 
halo profile is therefore given by 



^-3 ^=37TT '^'^ 

2 _ 2 / 2\ 

^^3 '^3e['-3) 

For e > 2/3, the density profile scales as 7--9«/(3<^+i)^ and is strictly stationary (s = 0). For e < 2/3, the 
density scales as r^^, and increases continuously with time. Moreover, we deduced also from the hydrostatic 
equilibrium equation that is the flattest stable density profile which can be obtained with purely radial 
orbits. 



3.2. Gas Halo Profiles 

We now turn to the pure gas collapse. The derivation that we presented in the previous section for dark 
matter can be applied to the gas dynamics, though the gas has an isotropic distribution in velocity space. 

Hydrostatic equilibrium writes for a gas fluid element M 

^ odP GM 

where P is now the pressure of the gas. Integrating again this equation between M and the currently 
shocked fluid element Ms{t) leads to 

P{M,t) = PS) + ^^M-, f x'-'^dx (19) 

where we use the same scaled variable as before x — M/Ms{t). Note that this last equation differs from the 
corresponding dark matter equation (eq. |l5|), due to the different dimensionality of the phase space. The 
integral converges or diverges now if p < 1/2 or p > 1/2. The case p < 1/2 has to be rejected in order to 
satisfy the hydrostatic equilibrium assumption. We then find the validity range for q 

e>l ,>0 

1 2/1 
^^6 '^^3^ 



The action takes the same form as for dark matter, and again, the least action principle applied on each gas 
fluid element leads to the self-similar solutions for gas 
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3e+ 1 
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7 
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3e 
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6 



1 



) 



For e > 1/6, the density profile scales as 7'-9f/(3's+i)^ and is strictly stationary (s = 0). For e < 1/6, 
the density scales as r~^, and increases continuously with time. These new solutions differ from the dark 
matter ones, although they show the same two characteristic regimes: a strictly stationary and hydrostatic 
flow for e > 1/6, and a non-stationary, but quasi-hydrostatic flow for e < 1/6. We also show that the flattest 
stable density profile is in that case. 



This analytical approach can be tested numerically. For dark matter, Fillmore & Goldreich (1984) 
solved the Vlasov-Poisson equations and described in a semi-analytical way the dark matter particles 
trajectories. They confirmed their analytical work, and therefore, the validity of the dark matter self-similar 
solutions. Bertschinger (1985) solved the pure gas collapse for e — 1, integrating semi-analytically the gas 
dynamics equations. He found that the gas density profile was similar to the corresponding dark matter 
profile. We recovered also this result in the previous section. 

To test our gas self-similar solutions, we use here the spherical hydrodynamical code, presented in 
Chieze, Teyssier & Alimi (1996), to which the reader is referred. Shock waves are followed using the 
pseudo-viscosity method, with a tensorial formulation of the viscous stress (Tscharnuter & Winkler 1979; 
Chieze et al. 1996). The shock front is captured within two or three cells, without post-shock oscillations. 
These calculations allow us to recover the asymptotic halo regime, and also the exact location of the shock 
front. The self-similar initial conditions are introduced in our code as follows: the general power-law 
density profile is smoothly connected to a non-singular, inner core, containing less than 0.1% of the final 
turn-around mass. We then choose an initial epoch, such as the density contrast in the center is equal to 
1%. The initial gas temperature is uniform and equal to the cosmic background radiation temperature. 
This is in any case a very small fraction of the final halo temperature (typically 10~^). The hydrodynamical 
equations are then solved using physical coordinates (r, t), rather than the self-similar coordinate X{t). 

In figure (|^), we plot the mass, density, velocity and pressure profiles obtained for the e = 1 run. We 
compare our results with the semi-analytical profiles obtained by Bertschinger (1985). Note that both 
results agree remarkably well. Moreover, the computed shock front position is recovered within 3% accuracy. 
In figure (H), we plot the density profiles obtained for e = 0.8, 0.6, 0.4 and 0.05. The straight line shows in 
each graph the corresponding asymptotic power law derived analytically in the last section. The numerical 
results agree perfectly well with the analytical ones. The e = 0.05 case illustrates the change of behavior in 
the halo dynamics (e < 1/6). The power law p oc l/r is very well recovered. We list in table (|l]) the shock 
front position computed by our hydrodynamical code. Note that for lower e, the shock front is located 
at a much deeper radius. 



4. 



NUMERICAL SOLUTIONS 
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5. DISCUSSION 

We now discuss possible applications of our results on cosmic structures formation theory. We discuss 
here three major consequences of our work : an origin for segregation between gas and dark matter in X-ray 
clusters (David et al. 1996), an explanation of the Navarro et al. (1996) density profile and a justification 
for the stable clustering hypothesis. 

5.1. Mass Distribution in the Core of X-Ray Clusters 

The first application of the present work is to propose an interpretation of Navarro et al. (1996) dark 
matter density profile found with high resolution N-body simulations (see eq. [0). This profile scales as 
in the center of the cluster, as for a large range of intermediate radii and as in the outer regions. 
The outer regime, located roughly at the Virial radius, is likely to be the relaxation layer where shell 
crossings dominate, which leads to a profile steeper than isothermal. The intermediate power law (r^^) is 
expected from the analytical approach for e < 2/3 and for purely radial orbits. Moreover, Tormen, Bouchet 
& White (1996) recently studied the distribution function in velocity space of simulated dark matter halos. 
They calculate for that purpose the factor f3{r) = 1 — CT|/2cr^ where ut and ct^ are respectively the tangential 
and the radial velocity dispersion at a given radius. /3 = 1 corresponds to purely radial orbits; /3 = to 
isotropic orbits and (3 — ~oo to purely circular orbits. In the region where the density scales as r^^, they 
found /3 ~ 1, while in the central region, where the density scales as r~^, they found /? ~ 0. In this paper, 
we show that for /? = 1, the flattest stable density profile is r~^, and for the gas, which has (3 — 0, the 
flattest stable density proflle is . The results presented in this paper suggest that the change of slope in 
the density proflle may be due to a change of dimensionality in the dark matter velocity space. 

5.2. Segregation between Gas and Dark Matter 

We showed that for e < 2/3, the density profiles obtained in the pure dark matter case differ from 
those obtained in the pure gas collapse. In cosmic structures, we expect both fiuids (gas and dark matter) 
to be dynamically important. Therefore, the coupled dynamics has to be studied in order to have definitive 
answers. We have not been able to solve here analytically the coupled case, mainly because self-similarity 
is broken by the presence of two distinct fluids. The total mass proflle is however expected to have a 
power law, between the pure gas and the pure dark matter proflles with 3/(1 + 3e) < 7 < 1. The exact 
value of 7 depends on the value of Ub, the universal gas fraction. In Chieze et al. (1996), we study the 
case of the single spherical Fourier mode. It can be shown that this corresponds to e = 4/9. As expected, 
we obtain a segregation which depends on the chosen value oi Hb. Bertschinger (1985) studied the case 
e = 1 in the coupled case, assuming a very low value for Qb- As expected, he found no segregation. Our 
belief is that there is always a segregation between gas and dark matter for e < 2/3, and no segregation for 
e > 2/3. This would be due to the fact that the dark matter evolution is not strictly stationnary in the 
halo. This results in a systematic delay of gas shells relative to dark matter shells. Following Hoffman & 
Shaham (1985), it is possible to flnd the most probable initial density proflle around a peak in the initial 
Gaussian random fleld. This leads to a power law initial perturbation which is roughly proportional to the 
linear two-point correlation function with e — (n + 3)/3. We propose here that power spectra with n < — 1 
lead to segregation between gas and dark matter, with a dark matter density proflle p oc and a gas 
density proflle slightly shallower. However, if n > —1, gas and dark matter are expected to have similar 
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density profiles with p oc 7'-9«/(i+3e) This rather strong prediction has to be tested in 3D hydrodynamical 
simulations. Anninos & Norman (1996) found no significant segregation in their 3D simulations. They used 
a standard CDM power spectrum where — 1 < n < on clusters scales. Moreover, they found that the halo 
density profiles are well fitted by a r~^/* power law, both for gas and dark matter. 

5.3. Stable Clustering Hypothesis 

The last consequence of our work concerns one of the most powerful analytical prediction in 
gravitational dynamics, namely the stable clustering ansatz (Peebles 1980). This hypothesis states that on 
small scales the relative pair velocity cancels exactly the Hubble flow. For scale-free power spectra, this 
leads to the non- linear two-point correlation hmction ^ cx r"'*' with 7 = (9 + 3n)/(5 + n). This regime is 
valid for ^ 1 and has been extensively tested by N-body simulations (Jain 1996). Padmanabhan et al. 
(1995) try to find another derivation of the stable clustering hypothesis using the secondary infall model 
(which is the main topic of our paper). Instead of the linear two-point correlation function, they used for 
the initial spherical perturbation the smoothed r.m.s. density contrast. This leads to e = (n-|-3)/6. Using 
the solution of Fillmore & Goldreich (1984), the resulting halo profiles are (Jain 1996) 

n > 1 p cx r-(9+3n)/(5+n) 

n < 1 p (X 

The next step is to assume that the two-point correlation function in the highly non-linear regime 
is dominated by the density field of such spherical halos. This has led Jain (1996) to conclude that the 
stable clustering hypothesis might be invalid for n < 1. However, Fillmore & Goldreich (1984) results apply 
only for purely radial orbits. Therefore, the stable clustering regime might be recovered within a fully 3D 
density field. We saw in section 5.1 that recent N-body experiments suggest that on very small scales the 
distribution of dark matter particles is nearly isotropic in velocity space. Therefore, to derive the halo 
profiles using the initial conditions of Padmanabhan et al. (1995), one has to use the new solutions we 
found in this paper. This reads 

n > -2 p oc r-(9+3")/(5+») 
n < —2 p oc 

Therefore, assuming that the two-point correlation function is dominated in the highly non-linear 
regime by the density profiles of quasi-spherical halos with an isotropic distribution of particles in velocity 
space, one could justify the stable clustering hypothesis for n > —2. Significant deviation could be detected 
for n = —2.5, but this rather fiat power spectra might be difficult to study using N-body simulations. 

5.4. Conclusions 

In this paper, we derive analytically gas and dark matter self-similar solutions for the gravitational 

collapse of a spherical perturbation embedded in an expanding Universe. The gas self-similar solutions can 
be easily extended to dark matter particles with an isotropic distribution in velocity space. The behavior of 
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the gas differs from the dark matter one, leading to different results for each fluid. We find three possible 
consequences of our work : a possible segregation between gas and dark matter for a certain class of initial 
conditions (e < 2/3), a possible explanation for a dark matter density profile p oc in the central region 
of cosmic structures, and a possible justification of the stable clustering hypothesis for power spectra with 
n > -2. 
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e 


1 


0.8 


0.6 


0.4 


0.2 


0.05 


A. 


0.36 


0.33 


0.29 


0.23 


0.13 


6.2x10-3 


AA, 


0.01 


0.01 


0.01 


0.01 


0.02 


3.0x10-3 



Table 1 : Self similar radius of the shock front and the corresponding uncertainties due to finite resolution 
obtained for different values of e by our hydrodynamical code. 
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X 




Fig. 1. — Density, pressure, velocity and mass profiles obtained in the e = 1 case for the pure gas collapse, as 

a fonction of the self-similar radius A = r/r*(t) (sec text). Solid lines are the results obtained by Bertschinger 
(1985) using a semi-analytical method. Open circle are the results obtained by our hydrodynamical code. 
Both results agree perfectly well. 



- 14 - 



10- - 



10' 



10^ - 



10" 



e=0.8 



-36/17 




10° 



10" 



Q 10* 



10' 



10-^ 



10-^ 



10' 



10" 



10" 




e=0.4 



10' 



10-^ 



10-^ 



10-' 



10" 



10"^ 
10" 
10^ 
Q 10' 
10* 
10^ 
10° 





1 

-1 




r 


-9/23 




r 




e=0.05 


1 



10-'" 10° 10° 10"' 10" 

X 



10" 



10° 



10° 



10' 



10^ 



10" 



10" 





1 

-9/8 

r 


e=0.2 


1 



10 10 10 10 10 10 10 



Fig. 2. — Density profiles obtained for e =0.8, 0.6, 0.4 and 0.05 by our hydrodynamical code, as a fonction 
of the self-similar radius A = r/r^{t). The straight line in each graph is the power law predicted by the 
analytical calculations performed in this paper. 



